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We investigate the persistent currents, spin-polarized local density of states, and spectral func¬ 
tions of topological superconductors constructed by placing ferromagnetic impurities on top of an 
s-wave superconductor with Rashba spin-orbit interaction. We solve self-consistently for the super¬ 
conducting order parameter and investigate both two-dimensional blocks and one-dimensional wires 
of ferromagnetic impurities, with the magnetic moments pointing both perpendicular and parallel 
to the surface. We find that the topologically protected edge states of ferromagnetic blocks give rise 
to spin-polarized edge currents, but that the total persistent current flows in opposite direction to 
what is expected from the dispersion relation of the edge states. We also show that the Majorana 
fermions at the end points of one-dimensional wires are spin-polarized, which can be directly related 
to the spin-polarization of the edge currents of two-dimensional blocks. Connections are also made 
to the physics of the Yu-Shiba-Rusinov states for zero-dimensional impurities. 
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I. INTRODUCTION 

During the past one and a half decade, topology has 
come to play an important role in several different but 
related aspects of condensed matter physics. In 2001, 
Kitaev proposed that a one-dimensional (ID) spinless p- 
wave superconductor can provide a route towards topo¬ 
logical quantum computation, since it can host topologi¬ 
cally protected Majorana fermions.® A few years later the 
subject of topological insulators was born, with the es¬ 
sential feature being an insulating bulk and topologically 
protected metallic interface states inside the bulk gap.®® 
The close connection between Majorana fermions and the 
edge states in topological insulators was explicitly made 
clear through the applic ation of topological band theory 
to superconductors .1121111 An early proposal of such a topo¬ 
logical superconductor involved interfacing a topological 
insulator with a conventional s-wave superconductor.!^ 
This was followed by proposals to construct a related 
setup using more conventional building blocks, such as 
s-wave superconductivity, magnetism, and Rashba spin- 
orbit interaction.E^HUlxhis latter approach is particularly 
interesting because it can physically be realized in many 
different ways and therefore leaves plenty of room for ex¬ 
perimental implementation. 

One of the most recent experimental investigations of 
Majorana fermions involves putting ferromagnetic impu¬ 
rities on top of a superconductor to create a ID topologi¬ 
cal superconductor wire, where surface effects give rise to 
the needed Rashba spin-orbit interaction.^^St^l] xhe same 
setup can also be used to assemble structures such as 
2D ferromagnetic islands or blocks, which have attracted 
attent ion b ecause of their topologically protected edge 
states.!22l2ni With this particular approach in mind, but 
also knowing that the system can be implemented in 
a variety of ways,!23H^ we here investigate both block 


and wire configurations of magnetic moments on top of 
a Rashba spin-orbit coupled superconductor. To accu¬ 
rately model these systems, we solve self-consistently for 
the superconducting order parameter. 

In particular, we investigate the persistent currents, 
local density of states (LDOS), and the spectral func¬ 
tions. We recently investigated the current patterns for 
a magnetic point impurity,!^ providing insight into the 
properties of the individual atoms in ferromagnetic do¬ 
mains and also connecting to the well-known physics of 
Yu-Shiba-Rusinov states of magnetic impurities in su¬ 
perconductors. Here, we instead start from the topo¬ 
logical band theory of a 2D ferromagnetic block and in¬ 
vestigate the resulting edge states, which are shown to 
be spin-polarized. We then study ID wires and connect 
their properties to the 2D blocks, by considering them as 
blocks with a width of a single site. Specifically, the end¬ 
point Majorana fermions are shown to be spin-polarized 
with a spin-polarization directly related to that of the 
edge states of the block. Through this approach we are 
able to fully elucidate the close connection between the 
topological states in 2D and ID, as well as the connection 
to the Yu-Shiba-Rusinov states in OD. In addition, we 
find persistent currents around the ferromagnetic block, 
which surprisingly flows in opposite direction to what is 
expected from the dispersion relation of the topologically 
protected edge states. This apparent contradiction is re¬ 
solved as the gap-crossing edge states are indeed found 
to give rise to a current compatible with their disper¬ 
sion relation, but lower energy quasiparticles are found 
to produce even stronger counter-propagating curren ts, a 
phenomenon known in other superconductors.!2^*29| This 
has important experimental consequences, since physical 
probes such as superconducting quantum interference de¬ 
vices (SQUIDs), which measure the magnetic field gener¬ 
ated from the total persistent currents, will give opposite 
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results to probes measuring transport properties involv¬ 
ing only states close to the Fermi level. 


II. MODEL 

Motivated by the presence of Rashba spin-orbit inter¬ 
action, s-wave superconductivity, and magnetic impuri¬ 
ties in recent experimental setups ,^^21111] consider the 
following Hamiltonian 

'H ='Hkin +'Hv^ +'Hso +'Hsc, (1) 

'Hkin = ~t ^ ~ M 'y ' 

(i.D.o- i.o- 

'^V^=-J2 • 0-)a<T' cLci<T', 

= I (e*®‘’4+b;CiT + H.c.) , 

ib 

T-Lsc = 5]] H.C.^ . 

i 

Here i and j are site indices on a square lattice, b is a 
vector pointing along the four types of nearest neighbor 
bonds, db is the angular coordinate of b, and cr is the 
vector of Pauli matrices. The parameters of the model 
are the strength of the hopping parameter t, chemical 
potential /r, Rashba spin-orbit interaction a, and su¬ 
perconducting pair potential 14c of the superconduct¬ 
ing substrate, as well as the effective Zeeman splitting 
coming from the ferromagnetic impurities ^ 2 ( 1 ). The su¬ 
perconducting pair potential is used to self-consistently 
determine the order parameter Aj = —KcicqCi^) = 
where and are the tr-spin 
electron and hole components, respectively, of the i^th 
eigenstate. We here set the hopping parameter t = 1, 
Rashba spin-orbit interaction a = 0.56, chemical po¬ 
tential /r = —4, and superconducting pair potential 
14c = 5.36. The specific values are not important, but 
our results are generally valid in a wide range of parame¬ 
ter space.® The Zeeman term is set to zero everywhere, 
except for a block- or line-shaped region, where it is set 
to a finite Vz, with the direction determined by the unit 
vector n. This allows us to study both 2D blocks and 
ID wires of ferromagnetic impurities embedded within a 
conventional superconductor. 


A. Currents 

The most essential property of topologically nontrivial 
systems is the appearance of gapless edge states. In the 
case of quantum Hall systems, these give rise to persis¬ 
tent currents along the edges, while for topological in¬ 
sulators they instead result in persistent spin-currents. 
Likewise, topological superconductors have gapless edge 
states, and it is therefore of interest to investigate the 


persistent currents in these systems. In fact, we find 
that a study of these currents help us understand also 
the spin-polarization of Majorana fermions at the end 
points of ID wires. We thus need not only expressions 
for calculating currents, but the spin-polarized currents 
for any given spin-polarization axis. 

To derive expressions for the currents, we consider the 
time rate of change of the spin-density operator pia- = 


( 2 ) 

Because we are interested in arbitrary spin directions, 
the operators have to be understood to be written in 
the basis of the particular spin direction of interest. If 
a is in the direction (0, p) in spherical coordinates, the 
operators are related to the operators in the basis of the 
Hamiltonian according to 

cL = <=0® (0 (0 

Now let ft denote the opposite spin direction to cr and 
define 


‘S'icr —Uicrcr Cjcr “t“ ^ 


Jibucr — ^Oib(TcrCj_|_j,^CiCT -b dibcro-Cig-Ci+bfr ) i 
? t I* t 

^ib(7(7 — “r ^ibacr^icr^i+bcr 


), (4) 


where ai^cr, Uibo-crj and Oibg-cr are coefficients to be deter¬ 
mined later. It is then possible to show, as is done in 
Appendix [b| that Eq. ([^ can be written as 


dpicT _ 

dt 



JihatT 


(5) 


where the sum runs over all bonds b for which there 
is a term in the Hamiltonian (for Eq. Q only nearest 
neighbor bonds). It is clear from Eq. Q that Sicr is an on¬ 
site source operator which converts d-spins into cr-spins. 
This operator is unrelated to currents and is thus of no 
further interest here. Next, Jibo-o- is a current operator for 
CT-spins moving away from site i along bond b, and is from 
Eq. ([^ seen to be responsible for decreasing the number 
of cr-spins on site i. Finally, also Jihaa is responsible for 
transferring cr-spins away from site i along bond b, but 
it flips the spin in the process. 

While Jibo-CT allows us to investigate spin-polarized cur¬ 
rents, the total current is obtained by adding all contribu¬ 
tions to the current and we are therefore also interested 
in the operator Ab = J2aa' 'dibo-o-'- Here the summation 
over ct' ensures that we add both spin-polarized currents 
as well as ’spin-flipping’ currents, while the summation 
over cr ensures that we consider the outflow of both types 
of spins. To be able to represent the current as a vector 
field, we also construct the on-site current operator by 
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adding the currents along the bonds away from site i into 
the directed-average vector operator Ji = 
and similarly for the spin-polarized and spin-flipping cur¬ 
rents. 

To fully construct the current operators we also need 
to evaluate the a-coefficients in Eq. Q. This is compli¬ 
cated by the fact that we are interested in spin-polarized 
currents for arbitrary spin-polarization axes. In addi¬ 
tion, the final expressions for the currents are obtained 
by evaluating expectation values (Ji) and (Jio-cr'): where 
the calculation also needs to be done in the particular ba¬ 
sis. A detailed derivation of the method used to calculate 
currents is provided in Appendix]^ In the following we 
present the current in two complementary ways. First of 
all, maxi |(Ji)| as a function of the magnetic moment al¬ 
lows us to identify the onset of persistent currents in the 
material. However, because the maximum is taken over 


the whole sample, it gives no information about where 
these currents appear. Therefore, for selected parame¬ 
ters we also show the current vector field. 


B. Spin-polarized LDOS 

For ID wires the topologically protected edge states 
manifest themselves in the form of Majorana fermions. 
For this reason we are interested in calculating the LDOS 
for a cut along the wire. As we want to relate the 
spin-polarization of these Majorana fermions to the spin- 
polarization of the edge currents of 2D blocks, we are 
further interested in the spin-polarized LDOS for a given 
spin-polarization axis. This is calculated using the ex¬ 
pression 






cos^d) 

cos(|) sin(|)e*‘^ sin^(|) 
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( 6 ) 


The total LDOS is obtained by adding the spin-polarized 
LDOS for two opposite spin-polarizations. 


III. DIMENSIONAL REDUCTION 

For a homogeneous bulk it is possible to write the 2D 
model of Eq. 0 aJi^ 


C. Spectral function 


where 


H(k) 


iJo(k) A(k) 

At(k) -ffo^(-k) 


( 8 ) 


The edge states of the 2D block also appear in the 
spectral function. Due to the limited size of the simu¬ 
lated system it is important to single out the contribu¬ 
tion to the spectral function coming from the edge states, 
otherwise bulk states pollute the edge spectrum. We do 
this by introducing a state classification function C{v), 
which is one for every state classified as an edge state, 
and zero otherwise. We here define a state to be an edge 
state if more than 50% of its probability amplitude is 
located inside a region, six sites wide and symmetrically 
positioned across the boundary. With this definition, the 
contribution to the spectral function coming from edge 
states can be calculated as 


A{\^,E)=Y,\u,k^^\'^5{E-E,)C{v), (7) 


where is the Fourier transform of u^xa, and u^xa 

is the restriction of Ui,i„ to the y-coordinate at which the 
bottom of the ferromagnetic block is located. 


i7o(k) =e(k) - Wa, - £o(k) • ct, (9) 

A(k) =iAay, (10) 

e(k) = — 2t [cos(/ca;) -I- cos(fcj,)] — /i, (11) 

£o(k) =a(sin{ky), - sin{kx),0). (12) 

This Hamiltonian have been topologically classified and 
the relevant c ondit ion for being in a topologically non¬ 
trivial phase iJSI^ 

0<\A\^ <V^-{4t + pf = V^, (13) 

where the last equality follows because p, = —At. More¬ 
over, making the restriction ky = 0 in Eq. 0 and defin¬ 
ing p = p + 2t, we obtain 

iLo'°(k) = - 2tcos{kx) - p - V^az + asin{kx)cry, (14) 

which describes a ID superconducting wire oriented 
along the x-direction, recently intensively studied in 
e.g. Refs. |17ll8l31j . It can easily be confirmed that this 
model, just like the 2D case, goes through a topological 
phase transition at |Ap = — {2t + pY■ This is one 
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example of a dimensional reduction from two to one di¬ 
mension, and similar relations play an important role in 
this work. 

Note, however, that if the ID wire is embedded in a 2D 
surface, this type of dimensional reduction makes little 
sense physically. First of all, setting ky = 0 corresponds 
to considering wave-functions that are completely delo¬ 
calized in the y-direction. Second, a ID infinite wire is 
produced by setting T4(i) = VzSy, not by ignoring the y- 
direction. A proper dimensional reduction would rather 
involve integrating out the y-dependence to arrive at an 
equation of similar form as Eq. (14|. This would poten¬ 
tially renormalize the parameters, and the condition for 
being in the nontrivial phase should then be understood 
in terms of these renormalized parameters. We explic¬ 
itly point this out here, because, while we have the strict 
condition in Eq. (13) to rely on for identifying the topo¬ 
logical phase transition of a 2D geometry, we have no 
such condition for a ID wire. Instead, we have to rely 
on a range of consequences, such as the appearance of 
Majorana fermions and other indicators related to the 
current, in order to identify the nontrivial phase of the 
wire. 


IV. RESULTS 

A. Block with magnetic moments perpendicular to 
surface 


We begin our investigation with a quadratic ferromag¬ 
netic block, with the magnetic moments perpendicular to 
the surface (n = z). In Fig. the maximum value of the 
current density and the size of the order parameter at the 
center of the ferromagnetic block is plotted as a function 
of the Zeeman term, which quantify the strength of the 
magnetic impurity moments. The condition for being in 
the topologically nontrivial phase, given in Eq. (131, is 
fulfilled in the shaded regions. It is clear that the maxi¬ 
mum current density is comparatively small in the trivial 
phase, raises up inside the topologically nontrivial phase, 
and then falls off as superconductivity is destroyed in the 
ferromagnetic block. We also plot the spatial distribution 
of the persistent currents in Figs. [2]|^for the three points 
marked 1-3 in Fig. For small Zeeman splitting the 
current flows anti-clockwise along the edge of the block, 
in agreement with a recent Ginzburg-Landau calculation 
for a circular disk geometry.l^ Next, in Fig. the per¬ 
sistent current in the nontrivial phase is seen to not only 
be much larger, but it also flows in the opposite direc¬ 
tion to that in the trivial phase. We also note that this 
current is located mainly inside the ferromagnetic region. 
Below we show that these currents are related to the ap¬ 
pearance of edge states, but surprisingly propagates in 
the direction opposite to what is expected from the edge 
state spectrum. However, before we move on to a more 
detailed study of the currents in the nontrivial phase, 
which is the main focus of this work, we also say a few 



FIG. 1: (Color online.) Bottom: Order parameter in the 
center of the ferromagnetic block (blue line) as a function 
of the Zeeman term 14 in the direction n = z. Absolute 
value of the Zeeman term is indicated with dashed lines. 
Shaded regions indicate the topologically nontrivial phase 
where |I4| > |A| > 0. Top: Maximum current density as 
a function of the Zeeman term. 



FIG. 2: (Color online.) Edge currents in the topologically 
trivial phase of a ferromagentic block magnetized along the 
direction h = z, as indicated by grey arrows. The current 
density is proportional to the arrow length, with the longest 
arrow corresponding to a current density of 2.8 • 10~'*. See 
also (1) in Fig. 


words about the currents in the third region. 

We see that the maximum current density remains 
fairly large even after superconductivity is destroyed in 
the ferromagnetic block region, and these currents are 
also located around the edge of the block. Strong edge 
currents therefore remain even after the superconducting 
gap has been destroyed, even though it is not possible to 
motivate their presence using topological arguments for 
the block region. A few spikes also appear in the maxi¬ 
mum current density, and these are due to the formation 
of vortex-like currents inside the block region, as illus¬ 
trated in Fig. 1^ The exact occurrence of these peaks as 
a function of the magnetization is dependent on the size 
of the block and these current patterns result from the 
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FIG. 3: (Color online.) Edge currents in the topologically 
nontrivial phase of a ferromagnetic block magnetized along 
the direction n = z. See also (2) in Fig.[^ 



FIG. 4: (Color online.) Currents at a current density peak in 
the collapsed gap state of a ferromagnetic block magnetized 
along the n = z. See also (3) in Fig. 

interplay of the Zeeman term, Rashba spin-orbit inter¬ 
action, and boundary conditions. A complete study of 
the currents in this regime is beyond the scope of this 
work, and we here only mention these effects to clarify 
the properties of the non-zero current in the third region. 


1. Reversed persistent current 

Having covered the basic behavior of the currents, we 
now investigate in detail the current in the nontrivial 
phase and relate it to the topologically protected edge 
states. In Fig. we plot the contribution to the spectral 
function from states classified as edge states for a cut 
along the lower edge of the ferromagnetic block. The left 
and the right figure corresponds to the nontrivial phase 
that occurs for negative and positive Zeeman term, re¬ 


spectively. It is clear that the edge states have opposite 
dispersion, in agreement with the two phases having op¬ 
posite Chern numbers. We now note that the slope of 
the edge state in panel (2) of Fig. implies that the 
edge state has a propagation direction to the right, since 
V ~ dE/dk. However, this is in contradiction to the di¬ 
rection of the persistent current plotted in Fig. which 
flows to the left along the lower edge. 

The apparent contradiction is resolved once the con¬ 
tribution to the current from the edge states is singled 
out. In Fig. we plot the contribution to the cur¬ 
rent from eigenstates with energy in the narrow interval 
—O.I < E < 0 around the Fermi level, and the current is 
seen to flow in opposite direction to the total current in 
Fig. [3 The somewhat unexpected conclusion is therefore 
that the persistent currents flows in opposite direction to 
the current carried by the topologically protected gap¬ 
crossing edge states. This have important experimental 
consequences as the measured current direction will de¬ 
pend on the physical probe used. If for example a SQUID 
is used to measure the magnetic field generated by the 
total persistent current, the current will be found to flow 
in opposite direction compared to the dispersion of the 
edge states. On the other hand, any transport measure¬ 
ment relying on excitations close to the Fermi level will 
find a propagation direction in agreement with the edge 
state dispersion. 

The reason for the reversed flow of the persistent cur¬ 
rents compared to the edge state currents can be un¬ 
derstood by considering why an ordinary bulk is free 
from persistent currents. This is not due to zero con¬ 
tribution to the current from individual charge carri¬ 
ers, but rather because the sum of all such contribu¬ 
tions cancel each other out. The total current car¬ 
ried by the system is obtained by integrating the prod¬ 
uct of the occupied charge carriers charge and velocity: 
J ~ /s(k)<fc^ eVkUdk = VkUdk = 0. Here 

the first step follows because the charge of each quasi¬ 
particle state is the same, while the second step follows 
because the spectrum is continuous and periodic in k. 
Neither of these assumptions are true at the edge of a 
topological superconductor. First of all, a branch in the 
spectral function which crosses the Fermi level an odd 
number of times cannot be described by a simultaneously 
continuous and periodic spectrum. Second, the quasi¬ 
particles in a superconductor are mixtures of electrons 
and holes and do not all carry the same charge. Consid¬ 
ering the persistent current carried by the edge states, 
we now note that only half of the branch corresponds 
to occupied states. Further, the occupied edge states 
are strong mixtures of electrons and holes because they 
are close to the Fermi level. The total current carried 
by this band is therefore much smaller than what would 
be expected if the band corresponded to a fully occu¬ 
pied electron band. It is this lack of occupied electron 
states associated with the gap-crossing band which fails 
to cancel the counter-propagating currents coming from 
the low energy states. This phenomenon of total currents 
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FIG. 5: (Color online.) Spectral function along the lower 
edge of a ferromagnetic block magnetized along the direction 
n = z. Only contributions from states classified as edge states 
are included. Edge state on the left and right corresponds to 
the topologically nontrivial phases at negative and positive 
Zeeman term, respectively, see (2’) and (2) in Fig. The 
opposite dispersions reflect the opposite signs of the Chern 
number in the two phases. 



FIG. 6: (Color online.) Same as in Fig. but when only 
contributions from eigenstates with energy in the interval 
—0.1 < E < 0 have been inclnded. 


cc-axis flow in opposite direction. However, in the edge 
regions the two currents flow in the same direction and, 
together with the ’spin-flipping’ currents, these gives rise 
to the total edge currents in Fig. 

We also note the partial depletion of a:-up spin cur¬ 
rents towards the right end of the block in Fig.j^a). The 
picture is radically different in Fig.j^b), where only con¬ 
tributions to the currents from states close to the Fermi 
level are included. Here, the a;-up spins are located on 
the right side, and with a similar argument to that in 
Sec. |IV A 1[ we can understand the depletion of the total 
cc-up spin current in Fig. [^a) as a partial, but incom¬ 
plete, cancellation by the states close to the Fermi level. 

It is clear from both the total current in the bulk re¬ 
gion and the edge currents in Fig. [^b) that x-up spins 
favor motion along the positive y-axis. This can be un¬ 
derstood as a consequence of the presence of Rashba spin- 
orbit interaction. The Rashba spin-orbit interaction cou¬ 
ples momentum and spin through k x cr, which leads 
to spin-momentum locking for the currents.l^ This spin- 
momentum locking and localization of different spins on 
different edges are, as we will discuss next, also important 
for understanding ID wire systems. 


B. Wire with magnetic moments perpendicular to 
surface 


We now move on to an investigation of a ID wire with 
the magnetic moments still perpendicular to the surface, 
which corresponds to the ID model in Eq. (14). As ex¬ 
plained in Sec. we no longer have a simple condi¬ 


tion, such as Eq. (131, to rely on to identify the topo¬ 


logical phase transition, because the system is embedded 
in a 2D surface rather than being truly ID. However, we 
find that this present little problem since we have several 
other means to identify the nontrivial phase. 


adding up from contributions of quasiparticle currents of 
states ’’inside the gap” and deep states that contribute 
to overall condensate flow is k nown in unconventional 
superconductors and superfluids! ^^ * ^^ ! 


2. Spin-polarized current 

Next we study the spin-polarization of the currents. In 
particular, the x and y spin-polarization axes are of in¬ 
terest because the Rashba spin-orbit interaction couples 
these axes differently to the momentum. In Fig. we plot 
the total spin-polarized current as well as that carried by 
eigenstates in the energy interval —0.1 < if < 0 for spins 
polarized along the x-axis. The spin-polarized currents 
for y-up, x-down, and y-down spins are identical to these 
under a spatial rotation of 7r/2, tt, and 37r/2, respectively. 
It is clear from Fig. [^a) that strong spin-polarized cur¬ 
rents flow through the block, but the total current in the 
bulk cancels, since the spins polarized along the negative 


1. Partial gap collapse and topological phase transition 

In Fig. we plot the maximum current density and the 
value of the order parameter at the midpoint of a wire. 
First of all, we note that compared to Fig.j^the horizon¬ 
tal axis has been rescaled and, compared to the block, 
a roughly twice as strong magnetization is required to 
suppress superconductivity on the wires sites. A similar 
calculation for a point impurity reveals that the magne¬ 
tization has to be about four times as large as as that 
of the block geometry in order to destroy superconduc¬ 
tivity. This can be understood from the fact that each 
edge site in a block is connected to a single nonmagnetic 
superconducting site, while each site in the wire and for 
the single impurity is connected to two and four such 
superconducting sites, respectively. 

Next, we note that the partial collapse of the order 
parameter in itself does not imply a topological phase 
transition, even though we numerically always see these 


























































7 



FIG. 7: (Color online.) Spin-polarized current for spins along the positive a:-axis for a ferromagnetic block magnetized along the 
direction n = z. (a): Total spin-polarized current, (b); only contributions from states with energy in the interval —0.1 < E < 0. 
The current pattern is related to the spin-polarized current for y-up, x-down, and j/-down spins by a spatial rotation of n/2, tt, 
and 37r/2, respectively. 



FIG. 8: (Color online.) Similar to Fig. but for a ID wire 
along the a;-axis and magnetic moments along the direction 
n = z. The topologically nontrivial phase is identified using 
the discontinuous onset of strong edge currents as well as the 
appearance of Majorana fermions. 


occurring at the same time. There are, however, good 
reasons for why these can be expected to occur at the 
same time. First of all, the topological phase transition 
occurs when A and 14 are related according to an expres¬ 
sion similar to Eq. (d), even though we can no longer 
rely on the last equality sign or the exact values of t and 
/i. The inequality is clearly valid when making the Zee- 
man term large enough, but also if the order parameter 
is small enough. A sudden drop in the order parameter 
will therefore likely induce the topological phase transi¬ 
tion. Thus, assuming no a priori knowledge of the value 
of A at which the transition occurs, the likelihood of the 
transition to occur is proportional to the size of the drop 
in the order parameter. 

An even stronger argument for why the two phenom¬ 
ena should occur in tandem can be made if we instead 
consider the drop in the order parameter to be induced by 


the topological phase transition. This relies on the fact 
that a topological phase transition is fundamentally re¬ 
lated to the interchange of energy levels across the Fermi 
level, also known as a band inversion. Because the order 
parameter is calculated using eigenstates below the Fermi 
level, the value of the order parameter can be expected 
to experience discontinuous changes each time a pair 
of eigenstates are interchanged across the Fermi level. 
Moreover, because the size of the order parameter in 
turn affects the eigenstates, this has enhancing feedback 
effects, which are only captured in fully self-consistent 
calculations such as those performed here. Discontin¬ 
uous changes in the order parameter and the topologi¬ 
cal phase transition are therefore intimately tied to each 
other. We note that this is closely related to the physics 
of single magnetic impurities. The associated Yu-Shiba- 
Rusinov states, which in non-self-consistent calculations 
cross the Fermi level continuously, cross discontinuously 
in self-consistent models because of the feedback effect 
of the occ upati on numbers of these states on the order 
parameter .E33I1 


2. Signatures of nontrivial topology 

So far we have discussed the relation between the col¬ 
lapse of the gap and the topological phase transition, 
arguing that there are strong reasons to expect them to 
occur at the same time. We have also mentioned that 
they always occur at the same time in our calculations, 
but so far not said what signatures we have used to iden¬ 
tify the topological phase transition. The first signature 
is the strong peak in the maximum current density in 
Fig.i Just like for the square block, the onset of the 
peak is also associated with a reversal of the direction 
of the persistent current, as can be seen in Fig. Most 
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FIG. 9: (Color online.) Persistent currents around ferromag¬ 
netic wire with magnetic moments along the direction h = z. 
Left: Trivial phase. Right: Nontrivial phase. See also (4) and 
(5) in Fig.[^ 


Infinite Finite Single site 



FIG. 11: (Color online.) Schematic figures showing the edge 
states along the rightmost edge of a ferromagnetic strip with 
varying length. If the edge is made infinitely long the spec¬ 
trum is continuous, but for any finite length, the spectrum 
becomes discrete. For a ID wire placed along the i-axis, the 
Majorana fermion at its end point can be understood as the 
single state that remains when the length of the i/-edge be¬ 
comes a single site. 


Qnaxi|(Ji)|=6^£^J^^^ 


^maxj |(ji)| = 5.8 ■ 10 



wire wire wire 

FIG. 10: (Color online.) LDOS as a function of position along 
a ID wire when magnetic moments are perpendicular to the 
surface (n = z). Both total and spin-polarized LDOS are 
shown, using the a:-axis as spin-polarization axis. The extent 
of the wire is indicated at the bottom. One Majorana fermion 
is clearly seen at each end of the wire. 


important though, this is also the point where Majorana 
fermions appear at the end points of the wire, as shown 
in Fig. 10 The Majorana fermion nature of these zero- 
energy states have been confirmed by checking that they 
appear as locally non-degenerate Bogoliubov-de Gennes 
quasiparticles. 


3. Majorana spin-polarization and block edge currents 

An interesting feature of the Majorana fermions seen 
in Fig. are their spin-polarizations. The Majorana 
fermions at the right and left end of the wire are polar¬ 
ized along the x-up and a;-down direction, respectively. 
That Majorana fermions have spin-polarization has al¬ 
ready been predicted,^ but we are here able to under¬ 
stand it in relation to the spin-polarized currents around 


a 2D ferromagnetic block. In Fig. Hb) it was shown 
that the contribution to the x-up polarized current com¬ 
ing from states close to the Fermi level is located on the 
right edge. This is in complete agreement with the x-up 
polarized Majorana fermion localized on the right end of 
the wire. This demonstrates the intimate connection be¬ 
tween the topologically protected edge states of a 2D and 
ID system. 

The connection between edge states and Majorana 
fermions can be made even clearer by considering the 
edge states plotted in Fig. For simplicity of argu¬ 
ment, we first consider an infinite ferromagnetic strip 
rather than a square block and let the infinite dimen¬ 
sion be along the y-direction. On each of the two edges 
one spectral branch cuts the Fermi level, one of which is 
schematically depicted in the leftmost picture of Fig. [m 
Next, consider the y-dimension wrapped around a cylin¬ 
der, such that it is still translationally invariant but with 
finite length. This results in a discretization of the spec¬ 
tral function as shown in the middle picture of Fig. [m 
By finally making the cylinder small enough such that the 
y-direction consists of a single site, a ID wire with a sin¬ 
gle Majorana fermion at each edge results, as depicted to 
the right in Fig. 11 This procedure closely resembles the 
dimensional reduction performed in Sec. |III[ and explic¬ 
itly shows the close relation between Majorana fermions 
in ID and edge states in 2D. However, this dimensional 


reduction does not, just like that performed in Sec. Ill 


correspond to the actual physical situation. This is im¬ 
portant to recognize because it might otherwise give the 
false impression that it is enough to discretize the edge 
spectrum to obtain a localized Majorana fermion. After 
all, the state at £1 = 0 in the middle panel of Fig. [D is 
also a Majorana fermion. 

The reason why the relation between the 2D block and 
the ID wire is slightly more complicated is that the two 
x-edges, which can be treated completely independent 
for the infinite strip, is joined by two y-edges. The spec¬ 
tral branches of the different edges therefore join into one 
single branch with eigenstates delocalized over the whole 
boundary. In particular, the Majorana fermions that oth¬ 
erwise would have existed at the right and left edge over- 
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lap and hybridize. It is only once the y-dimension be¬ 
comes small enough, such that the two y-edges start hy¬ 
bridizing and thereby gaps out the spectrum along these 
edges, that the two x-edges becomes independent of each 
other. It is therefore clear why it is not enough to have 
a discrete spectrum, such as for a block with finite edge 
length, to obtain localized Majorana fermions, and why 
a (quasi-)ID wire is required in spite of the close connec¬ 
tion between the two systems. 


C. Block with magnetic moments parallel to 
surface 



FIG. 12: (Color online.) Same as in Fig. but for magnetic 
moments along the direction h = x. 


Having investigated magnetic moments perpendicular 
to the surface, we now turn to a ferromagnetic block with 
the magnetic moments pointing along the x-direction 
(n = x). This system does not have a topologically non¬ 
trivial phase, but is of interest because of its relation to 
the corresponding ID wires. 

In Fig. we plot the order parameter and maximum 
current density. Interestingly, a strong peak is observed 
also in this case, even though there is no topological phase 
transition. However, the build up of the peak is not dis¬ 
continuous as in the previous cases. Instead, the maxi¬ 
mum current density is sizable already before supercon¬ 
ductivity starts to collapse in the ferromagnetic block. 
In fact, it is not only the maximum current density that 
is large, but also the total persistent current. This is 
clear from Fig.[^ where it is seen that the current flows 
through the whole sample, rather than being limited to 
the boundary of the ferromagnetic block. A similar be¬ 
havior is predicted using a Ginzburg-Landau model in 
Appendix 

The current flowing through the bulk region of the fer¬ 
romagnetic block can be understood as a consequence 
of x-up and x-down spins moving in opposite direction 
because of spin-orbit interaction, similar to what was 
seen in Fig. The addition of a Zeeman term in the 
x-direction creates an imbalance between the two spin- 
species and thereby also a total current through the sam¬ 
ple. The current continues to increase until the ferro¬ 
magnetic region enters the normal state. The peak in 
the total current density is therefore not due to a sudden 
onset of new physics, but rather because the continuous 
build up of current comes to an end as superconductivity 
is destroyed. We also note that the region following the 
partial collapse of the order parameter and up is com¬ 
plicated by the presence of vortex like currents inside 
the ferromagnetic bulk region. See for example Fig. 

For strong enough Zeeman term the currents becomes lo¬ 
calized on the left and right edges, similarly to what is 
predicted using a Ginzburg-Landau model in Appendix 

El 



FIG. 13: (Golor online.) Persistent currents before the partial 
collapse of the superconducting gap for a ferromagnetic block 


magnetized along the direction n = x. See also (6) in Fig. 12 



FIG. 14: (Color online.) Persistent current pattern at maxi¬ 
mum current density peak at (7) in Fig. 12 Vortex-like solu¬ 


tions continue to appear in the ferromagnetic region also for 
somewhat larger Zeeman terms, until the current eventually 
becomes localized on the left and right edges. 
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FIG. 15: (Color online.) Same as in Fig. 10 but for ferrom- 
ganetic wire magnetized along the wire (h = x and wire along 
x). Note that the Majorana fermions are spin-polarized along 
the 2 -axis rather than the x-axis. 


D. Wire with magnetic moments along wire 


We now consider a wire with the magnetic moments 
along the wire. This is the dimensional reduction of the 
block in the previous section, where the width of the 
block has been set to one along the y-direction. This 
system can, just like the wire in Sec. |IVB[ be described 
by a Hamiltonian of the form in Eq. (14), if the sub¬ 
stitution Gz —f Gx is made. This substitution presents 
no problem because the two Pauli matrices in Eq. (14) 
remain orthogonal, and the topological properties there¬ 
fore remain the same. It is in fact appropriate to think 
of the substitution as a rotation of spin indices around 
the 2 /-axis; Gx,Gz —>■ —In Eig. we plot the 
LDOS and spin-polarized LDOS using the z-axis as spin- 
polarization axis. It is clear that the Majorana fermions 
spin-polarization that previously showed up in the x- 
component is now present in the z-component. More¬ 
over, direct comparison to Fig. 10 reveals that the a:-up 
component goes into the 2 -down components and vice 
versa, exactly as implied by the rotation of spin-indices. 


Finally, the maximum current density and the value of 
the order parameter at the wires midpoint are shown in 
Fig. |16| The strong discontinuity in the maximum cur¬ 
rent density is correlated with the appearance of Majo¬ 
rana fermions. We also mention that the current pattern 
for the wire closely resembles that for the ferromagnetic 
block in Fig. The current flows through the wire 
perpendicular to the direction of the wire, and the two 
patterns are so similar that we refrain from showing it in 
a figure of its own. 



FIG. 16: (Golor online.) Same as in Fig. [B but for a ferro¬ 
magnetic wire magnetized along the wire (n = x and wire 
along x). The topologically nontrivial region has been identi¬ 
fied using the discontinuous onset of strong edge currents, as 
well as the appearance of Majorana fermions. 


xlO 



FIG. 17: (Golor online.) Same as in Fig. but for wire with 
in-plane magnetization perpendicular to wire (n = x and wire 
along y). 


E. Wire with in plane magnetic moments 
perpendicular to wire 


The final setup we consider is a wire with in plane mag¬ 
netic moments perpendicular to the wire. Similarly as 
above, this system can be obtained by making the width 
of the ferromagnetic block in Sec. |IV C| equal to one in 
the cc-direction. However, because the dimension is re¬ 
duced in a different direction, the substitution of Pauli 
matrices needed to obtain an effective equation of the 
form in Eq. is different. This time the substitution 
kx,Gy,Gz —>■ ky,—Gx,Gx is required, and the two Pauli 
matrices for the spin-orbit and Zeeman term become the 
same. The consequence is that the system has no topo¬ 
logically nontrivial phase, and no Majorana fermions ap¬ 
pear. The maximum current density and order parameter 


at the midpoint of the wire is presented in Fig. 17 
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V. SUMMARY 

We have studied ferromagnetic impurities on top of a 
conventional s-wave superconductor with Rashba spin- 
orbit interaction. We find a close relation between the 
topological properties of 2D blocks and ID wires of the 
ferromagnetic impurities, as well as the physics of Yu- 
Shiba-Rusinov states of OD point impurities.!^ 

When the magnetic moments are directed perpendic¬ 
ular to the surface, both the 2D block and ID wire has 
topologically nontrivial phases. For a 2D block the non¬ 
trivial phase gives rise to spin-polarized edge currents. 
Interestingly, the total persistent current flows in op¬ 
posite direction to the current contribution from the 
topologically protected gap-crossing edge states. This 
means that different experimental probes will observe 
different current directions, depending on if they mea¬ 
sure the total persistent current or only that carried by 
low-energy excitations. The Majorana fermions that ap¬ 
pear at the end points of ID wires are spin-polarized, 
and their spin-polarization direction are directly related 
to the spin-polarized edge currents of the 2D block. The 
spin-polarization of the Majorana fermions is along the 
direction of the wire. 


When the magnetic moments are parallel to the sur¬ 
face the ferromagnetic block does not have a topologi¬ 
cally nontrivial phase. However, ID wires can still have 
a topologically nontrivial phase, depending on if the wire 
is aligned with the magnetic moments or not. When the 
wire is parallel with the magnetic moments, a topolog¬ 
ically nontrivial phase exist and the corresponding Ma¬ 
jorana fermions are spin-polarized along the z-axis. For 
both the ferromagnetic block and the two wire directions, 
a substantial amount of current flows through the bulk, 
as long as superconductivity is present in the ferromag¬ 
netic region. 
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Appendix A: Ginzburg Landau approach 


To complement the self-consistent calculations on a 
lattice, we here also provide the results of a Ginzburg- 
Landau (GL) treatment. The GL method gives an ap¬ 
proximate description of the system, but is unable to cap¬ 
ture microscopic details such as spin-polarization, edge 
states, density of states etc.. Nevertheless, it may be 
interesting to compare the GL results with the exact nu¬ 
merical diagonalization discussed throughout the paper. 
We consider the following dimensionless GL free energy 
functional 


F= 


\ (Al) 


-tfclVV'P + - -{zxS) 


Here, ipi'f') describes the superconducting order parame¬ 
ter in a 2D box, i.e. r = {x,y) and —L/2 < x,y < L/2, 
where L = 1. The magnetic block of size d = 0.3L is 
placed in the center with the magnetization along the x 
axis, i.e. S{r) = Sx for —d/2 < x,y < d/2. The fer¬ 
romagnetic vector S is analogous to the Zeeman term 
Vz used throughout the main text. The first term in 
Eq. (Al) is conventional and favors the superconduct¬ 
ing order at \ip\ = 1, whereas the second term describes 
the diamagnetic Pauli energy disfavoring superconduct¬ 
ing order. The superconducting order is destroyed when 
the Pauli energy exceeds the the superconducting con¬ 
densation energy (the so-called Chandrasekhar-Clogston 
limit), which occurs for S > Sp = 1 in Eq. (Al). In 
order to smoothly interpolate the Pauli term between 
the superconducting and metallic phases we introduced 
a sigmoid function 


as 


dtp 6F 

dr dtp 


(A3) 


Using Eq. (|A3|), the order parameter is updated itera- 

-'S-FUo 


tively as V^new = l/'oid - dx- 


Sip 


The time evolution 


stops for the solution tp{r) satisfying = 0, which, 

thus, minimizes the free energy. 

In Fig. 18 we plot the order parameter magnitude \tp\, 
phase Arg^l, and the current 


j = ib {tP*Vtp - VtP*tP) + c\tp\^ {zxS). (A4) 


The panels in the top (bottom) row correspond to the 
field S = 0.9 (^ = 2), which is below (above) the Pauli 
limit Sp = 1. For the smaller field (top row), the super¬ 
conducting order \tp\ is uniform, whereas the supercon¬ 
ducting phase is non-uniform and is distributed in such a 
way that the total current is continuous. In contrast, 
for the large field (bottom row), the superconducting 
order is destroyed in the region covered with the ferro¬ 
magnetic block. The superconducting phase is strongly 
uniform outside the ferromagnetic block. The current 
is suppressed everywhere, except in small vicinity of the 
ferromagnetic block. 

In general, the GL modeling agrees qualitatively with 
those discussed in Section |IV G| the current grows lin¬ 
early with the magnetic field until superconductivity is 
destroyed. At the same time, we are not able to capture 
the non-linear peak of the current shown in Fig. |12| as 
well as the appearance of the vortex structure in Fig. |14| 


Appendix B: Derivation current expressions 
1. Introdnction 




V = b. 


(A2) 


Note that ~ 1 (0) in the absence (presence) of the 
superconductivity, i.e. if tp rv 0 {tp ^ 0). To avoid pos¬ 
sible spurious numerical effects, we also add a boundary 
term with a large coefficient a{r) ~ 100, which destroys 
superconductivity at the boundary of the 2D box. The 
terms on the second line of Eq. (Al) describe the kinetic 


energy. The term proportional to & = 0.0001 describes 
the conventional superfluid kinetic energy, which defines 
the characteristic coherence length ^ ^ y/b — 0.01. The 
last term, proportional to c = 0.003, describes the mag¬ 
netoelectric couplin^^ between the current and magne¬ 
tization S. 


We search for the minimum of the free energy (Al) 


numerically. We employ the gradient descent method in 
which the order parameter evolves in imaginary time t 


Here we derive the general expressions for calculating 
the currents in an arbitrary spin-basis. The imple¬ 
mentation of the numerical calculations are based on 
the algorithm described in Sec. |B 6| The method itself 
is fairly straightforward. However, the derivation is 
complicated by the fact that the spin-polarized currents 
are not necessarily calculated in the same basis as the 
Hamiltonian. A short summary of the appendix follows: 

\B Notation 

Introduction of useful notation and derivations of 
important identities. 

\B 3[ ' Generalized current operators 

The main problem is formulated. General sink-source 
and current operators are defined. In summary, there 
are three important terms for any given spin: The sink- 
source term Sg , the current on bond b, and the spin- 
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FIG. 18: (Color online.)(a, d) Order parameter magnitude \ip\, (b, e) phase Argip, and currents calculated using Ginzburg- 
Landau theory. The position of the ferromagnetic block is indicated with a red dashed rectangle. The ferromagnetic vector is 
in-plane S = Sx. Top row corresponds to 5 = 0.9 and bottom row to 5 = 2. The magnitude of the currents in panel (f) is 
six-fold magnified compared to panel (c) as indicated by the magnification ratio shown in these panels. 


flipping current j^g. on b, which together satisfies 

^ = Ss-J2 ■ (Bi) 

b 

To derive these operators, we need to find a certain set 
of coefficients 


and Rashba spin-orbit interaction. The total table for 
a Hamiltonian containing all these terms is also provided. 

\B Note on superconductivity 

It is shown that the introduction of an s-wave super¬ 
conducting term in the Hamiltonian does not modify 
the method described in the rest of this appendix. 


\B 4\ Arbitrary basis 

General expressions are derived for how to calculate 
the g coefficients in an arbitrary basis from the 

knowledge of those same coefficients in the basis of the 
Hamiltonian. 

\B 51 Calculate sink-sources and currents 

Expressions for the expectation values of the sink- 
source and current operators are derived. 


2. Notation 


Here we introduce useful notation and show relations 
which will be important later on. This section is purely 
algebraic, and no physical motivation will be provided. 
The usefulness of the expressions will only become clear 


once general currents are treated in Section B 3 


\B Summary of method 

An algorithmic summary of the method is provided. 


a. Bond decomposition of commutators 


\B 7| ' Tabulating relevant a’s 

The exact number of g that needs to be 

evaluated in the basis of the Hamiltonian are deter¬ 
mined, and a useful way of tabulating them is introduced. 

\B <$[ • Evaluating tables 

Tables are evaluated for the following types of terms 
in the Hamiltonian: chemical potential, kinetic, Zeeman, 


Consider a general quadratic Hamiltonian 

H = Y.h%, 

otXX' 


where 


(B2) 


(B3) 
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Here i,j are site indices, A, A' are spin indices, and a etc., 
is an index distinguishing between different terms in the 

Hamitlonian such as kinetic, chemical potential, Zeeman, We now consider expressions of the form 


i 

h 


[^AA') 


(^^js5\'Kcl^CsK' SisS\K,icli^Cjy^ 
ij 

y ^ /*sA,s+fe,A'*^A k'. 


(B4) 


The summation in the last expression is said to run over 
all bonds b connected to site s, where a bond between site 
s and s + 6 is said to exist if the Hamiltonian contains a 
product of creation and annihilation operators with these 
two indices. In particular, the set of bonds has to be 


understood to include the bond between a site and itself if 
any on-site terms are included in the Hamiltonian. Using 
Eq. (B4), it is possible to decompose the commutator for 
the complete Hamiltonian as 


i 

h 




i 

h 


[^AA'>cLcg„/] 

c^AA' 


i 

h 


E 




Defining 


XX' kk' 


— Sx-q^K'fl'bx'K,^ 'Y, /z,A,7,A'i 


^XX' kk' 


— bK.r]b\'rj'5xK' ^ 


the expression can further be written as 


(B6) 


(B5) 


i 

h 




= E 

bXX'r}r}' 


XX'kk' t „ 
^s+b,rj,s,ri'^s+b,ri^^ 


_ lAA'kk' U ^ . A 

^s,ri,s-\-b,r)' ^sr)^s-\-b,r]' j 


-n 

bX 


Xkkk' t „ 
“s-|-b,A,s,K'S+b,A^s 


, _ bKAftfc' t r , 1 , 
SK '^S,K,S + b,X''SK''S+b,X 


(B7) 


b. Useful expressions 

The Hamiltonian is Hermitian if B = W , which using 
Eq. ( |B2[ ) implies 

E f^%x'4xCjy = Y (B8) 

OiXX'ij aXX'ij 


Identification of coefficients leads to t he co nclusion that 
Ea /a,A' = Ea ffxHX^ and using Eq. (@ we have 


uXX'KK X' Xk' K* 


(B9) 


It is therefore possible to rewrite Eq. (B7) as 
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i 

h 


[-^jCg^CsK'] — ^ ^ 

b\\'7)rj' 


f W'kk' 
\^s-\-b,ri,s,r}^ ^ 


s+b,r}^s7]' 


+ CL 


X'Xk'K* t ^ 

s-\-b,r)' ,s,r]^sr)^s-\-b,r) 


') 


^i:( 

bX 


Xkkk' t „ 
^s+&,A,s,k' 


Ak'k'k* t 
^s+&,A,s,k^sa 


,^s+b 


.a), 


(BIO) 


which eliminates the need to calculate b coefficients in¬ 


dependently. We also note from Eq. ( |B6| ) that 


XXkk, 


which together with Eq. (B9) gives 


AAKK _ Q 


(Bll) 


Using k! and r]' to denote the spins which are opposite to 
k' and ry', respectively, we also find from Eq. (B6) that 


XX' kk' _ XX' kk' 


c. Obtaing a’s through projection 


(B12) 


In Eq. (B6) explicit expressions for how to calculate 


were provided. It will also be useful to have a 
way of extracting coefficients from a general quadratic 
expression, also when the explicit expression for these 
coefficients are not known. To this end we consider a 
general quadratic expression 


= E- 


u, 


(B13) 


where fi and v are some arbitrary set of indices. It is 
clear that 


Cpecj),|0) = Qpp'. 


(B14) 


Applying this to Eq. (BIO) we can for b ^ 0 write 


Cs^b,7] J. Cgj^/|0) — ^ ^ ^5+6,?7,s,y?^ “ 

AA' 


a 


5 + 6 , 77 , 8 , 77 ' ’ 

(B15) 


while for 6 = 0 we have 


+ al-;\ (B16) 


which can also be written as 


k [+Cs^Csk;'] cJ^/|0) — 


a 

TjKKK,' 


0 

7]K,K,k' 

^sr}sr}' 

7)' k’ K.' K,*- 
ST)' ST] 


\i T] ^ n and t]' ^ k,' , 
\i T] ^ n and r]' = n' ^ 
\i T] n and v( ^ ^ 


,. I 7 ] K K. K*^ ‘ r J/ / 

^s,n,s,r,' + if ?7 = K and ry = k. 


r 


(B17) 


3. Generalized current operators 


The right hand side of this expression is of the form con¬ 
sidered in the previous section, and we therefore know 


With the notation introduced in the previous section, 
we are now ready to formulate the problem in a way 
that eventually will allow us to calculate spin-polarized 
currents in any given basis. To this end we consider the 
local density operator for particles of spin a on site s, 
which can be written as psa = clo-Cso-- The time rate of 
change of this operator is 


t^Psa 

dt 


i 

h 




(B18) 
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from Eq. (BIO) that it can be written as 


These can all be written as 


dpsa 

dt 




bX 

\ (J S (J ' rr Q rr ! 


-*'8,(7,8,(7/ 

acrcro- J c- 

^S,a,S,a^SG^SG ' ^ rr^<irr^SG 


^s,g,s,g'^sg'^sg 


+ (^^s+b,a,s,a'^l+b,a-^S(j + 0,s+b,a,s,adla^s+b,(7^ 

MO 

+ g.Cscr + fls+6,CT,s,(TCso.Cs+6,CT^ ) 


(B19) 


where it is the opposite spin of a. In the last step the 
on-site terms have been singled out of the expression, 
and the remaining sum over the bonds have been divided 
into one part free of spin-flips, and one involving spin- 
flips. Using Eq. (Bll) we can further conclude that the 
first term is zero. The next two terms converts a- and 
tt-spins into each other on-site. When considering the 
density of a single spin species such as a, these terms 
therefore acts as sink and source terms, respectively. We 
therefore define the sink-source operator on site s as 


C _„o-o-cr(T t „ 

—'^s,a,s,a'^sa'^ea 


8,(7,8,(7^8(7^ 


(B20) 


The third line describes processes where a cr spin jumps 
between site s and site s -I- 6. It is therefore the cur¬ 
rent into site s along bond b. We accordingly define the 
current operator along bond b as 


(„<y<yry<T t 

\^^8 + 6,(7,8,(7^S + 6,(7 


Cs(7 


I (7(7(7(7* n I 

“1“ ^8 + &,Cr,8,7-^87-^S + 6,7- J 7 


(B2I) 


where it is implicitly understood that the operator is only 
defined for actual bonds between different sites. We have 
here included an additional minus sign in order to define 
the direction of the current to be along the bond away 
from site s. The third term is similarly related to the in 
and out flow of cr-spins on site s. However, in contrast to 
the ordinary currents the spins flip during the transition. 
We therefore define the spin-flipping current along bond 
b as 


Jl 


b 

GG 


( aaaa t 
y^S-\-b,G,S,G^S + b,G 


^SG 


I GGGG* 

' ^S-\-b,G,S,G 


Cs-\-b,c 


(B22) 


With these definitions we can now write 


dps<7 

dt 


b^O 



(B23) 


(B25) 

With knowledge of the eigenstates of the Hamiltonian, 
the expressions just derived allow us to calculate on-site 
sink-sources, as well as currents on bonds, using Sg = 
{Sg) and respectively. However, it is also 

useful to know the current on a specific site, rather than 
on the bonds. If b is the vector along bond b, from s to 
s -I- 5, we also define the two types of currents at site s as 

= (B26) 

^ b^O 

This is a directed average of all currents flowing out from 
site s and forms a vector field, in contrast to the scalar 
currents defined on the bonds. 


4. Arbitrary basis 


Although we in principle can proceed to calculate the 
relevant coefficients in Eq. (B24) using Eq. (B6), this 


requires the Hamiltonian to be expressed in the same 
basis as a. It would be a nuisance if each time we want 
to calculate the current in a new basis, the Hamiltonian 
would need to be rewritten in that basis, commutators be 
reevaluated, and finally also the expectation values (Sg) 
and (Mo-')- much better if the current in one basis 

can be evaluated from the commutators and expectation 
values calculated in the basis of the original Hamiltonian. 

Consider therefore operators with arbitrary spin direc¬ 
tion 


Cs(T —T d^Cg^^ 

Cga =dlcg^ - b*Cgi = b^Cg^ + biCgi, (B27) 

where a is some arbitrary spin, and a is the opposite 
spin. Further, f and j, here refers to up and down spins 
in the basis of the Hamiltonian. If a is along the direction 
(0, (fi), in ordinary polar coordinates relative to the basis 
of the Hamiltonian, the coefficients can be chosen to be 

=-^* =cos(^), 

6^=6; = cos(^)e-*^. (B28) 

and the relevant commutator can be written as 

^ [H, cLc..] = E Kd.' \ [H, cLc«.<] . (B29) 

kk' 


It is clear that in order to have a complete understand¬ 
ing of the currents in a given basis a, it is sufficient to 
evaluate the three coefficients 


J GGGG GGGG GGGG 

1.^8,O',8,(75 ^8 + fc,(7,8,(7 5 ^S-\-b,G,S,G 


(B24) 


It is clear that the left hand side can be calculated for any 
spin direction a if the right hand side has been evaluated 
in the basis of the Hamiltonian. However, we remem¬ 
ber that it is the expansion of the left hand side into a 
coefficients that allows us to calculate currents. Assume 
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therefore that the a coefficients in the expansion of the 
commutators on the right hand side have been evaluated. 
We denote these coefficients by a capital A, to distinguish 
them from the a coefficients in the expansion in the basis 
of the expression on the left hand side. 

We will now show how a coefficients in any basis can 


be calculated from a knowledge of the A coefficients. To 
this end, we remind ourselves that we in fact are only 
interested in calculating the three a’s in Eq. (B24|. Fur¬ 


ther, using Eq. (B15l and (B17), these can be written 
as 


={o\cs,,^ci\o) 


KK K K 


ZZ,s,. =(0|c«+,,.^cL|0) 


dt 


— (0|Cs-|-b,K" ^ [77, C^^Csk'] Cg^,„ |0) 

kk' k" k'" 

fTM....=(0|c«+b.a^cL|0) 

= 51 [77,cl^Cs«/] c^^,„|0). 


(B30) 


The two last expressions can be simplified straightfor¬ 
wardly by using Eq. (B15) and then summing over indices 


appearing in (5-functions in Eq. (B6), which gives 


^s-\-b 


^s-\-b,a 


■G _ T*|7 |2t \k" KKK 

^ \ ■^S-\-b,K" ,S,k' 

kk' k" 

,S,.= E (B31) 


5. Sink-sources and currents 


The on-site term is slightly more complicated due to the 
special cases in Eq. (B17|: 


E Kb^'b^"bl,„ together with Eq. ( [bM| ) and 

. ' ' L ' ^ / I—I ^ I i TT -r-i rf-VTTT t'TTVt 4~ ^ 4-1-*^^ yAT7-'Y-» cx 4" T ^*1^ T rcx 1 1 T 4-^*r» 


KK' K" K" 


Eq. (B321, we can now write the expectation values for 


(B32) the sink-source and current operators as 


(5.) + , 

kk' 

{Ji ,)=- E Kb^' +<?r.;..(cLc.+b,«')//), 

kk' 


(B33) 


Here also the expectation values on the right hand side which is indicated by the subscript 77. Using that 
have been expanded in the basis of the Hamiltonian, (4„Cs+b,„') = (c^^^,Csk)*, this can finally be written 
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=2Re J , 

\ kk' / 

{Jaa) — ~ 21^6 f fls+b,CT,s,cr ^k^k' {^s+b,K^‘>i^') H j j 

\ KK,' / 

{J<7s) = ~ 2Re I a^_|_f, g. (Cg_|_j, ,,CsK')_f/j • 


d. Calculate relevant a’s in the chosen basis 


(B34) 


s,a,s,a / ^ 

kk' k" k" 


h*h ,h nh* ( I 

Of^0^f0^/'0^,N y^s,K" ,s,k''' ^ -^s,k'" ,s,k" J 

^(7(T(T(T _ \ ^ u*\h Ak" KKk' 

^S + &,(T,S,(T ~ / ^ ■^S-\-b,K" ,S,k'•> 

KK' K." 

<ii%,s,a= E (B40) 


6. Summary of method 


e. Calculate sink-sources and currents 


Here we give a short algorithmic summary of the 
method. We note that the two first steps only needs to 
be performed once in the basis of the Hamiltonian. Re¬ 
peated application of the subsequent steps then allows 
for the sink-sources and currents in arbitrary bases to be 
calculated. 



a. Calculate expansion coefficients of commutators in the 
basis of the Hamiltonian 


(4^) = - 2Re 


^s-|-6,cr,s,(T / ^ 
kk' 




(B41) 


A XX' kk' 
■^ 5 + 6 , 77 , 8 , 77 ' 


^Xri^X' k^k't}' 


i 

h 


E fs+b,\,s,\>- 
a. 


(B35) 


See Sec. |B8| for already evaluated terms for chemical po¬ 
tential, kinetic, Zeeman, and Rashba spin-orbit interac¬ 
tion. 


b. Calculate expectation values in the basis of the 
Hamiltonian 


f. Current vector field 

From the currents on the bonds, a current vector field, 
defined on site s, is finally calculated as 

J+ = ^Eb+-') (B42) 

^ 6#0 

7. Tabulating relevant a’s 


^^s+b,K.^>>4 H ■ 


(B36) 


c. Choose spin direction 

Choose appropriate coefficients and b\^ for the spin 
direction cr of interest 

Csa = b^Cs\ + b^Csi- (B37) 

For a spin pointing in the direction (0, (/?), use 

b^ =cos(-), 

=sin(^)e"*‘^. (B38) 

The coefficients for the opposite spin follows as 

h =bl, 

h=-b*r (B39) 


Having summarized the method, it only remains to 
evaluate the a coefficients in the basis of the Hamilto¬ 
nian, which are denoted by a capital A. By counting the 
number of spin indices in Eq. (B6) to be six, it may seem 
like we in general need to evaluate 2® = 64 a’s per bond 
to be fully able to expand the commutators on the form 
in Eq. (BIO). However, the three (5-functions ensure that 
only 2® = 8 of these are non-zero. For any given bond, 
we can therefore restrict ourselves to evaluate s «' ■ 

Further, Eq. (B12) allows us to reduce this number to 
four, by fixing k' =t. The right hand side of Eq. (BIO) 
can then be completely determined for any site s from a 
full knowledge of the table 


(bj) 

, 4 nn 

Atm 

mn 

m.n 

m.n 

(1,1) 





(1,2) 










(n,n) 
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Here i, j are site indices, with n the number of sites. The 
elements necessary to reconstruct Eq. (BIO) for any given 
site s can be found in the table by setting i = s + b and 
j = s. Note that also most terms in this table will be 
zero. From Eq. (B6) it is clear that an a will only be non¬ 
zero if the Hamiltonian contains some f^xjX' connecting 
the two sites. It is therefore sufficient to further limit the 
evaluation of the table to those pairs of sites that share 
a bond. 

We also note that if we assume that each = 
Saa' ^AA' is Hermitian, all results derived in this section 
applies to the individual /i“’s. In particular, the entries 
in the table for the full Hamiltonain H can be obtained 
by adding the corresponding entries in the partial tables 
of the individual /i“’s. 


where 


flx,x'=5xx'Ur (B46) 

The table is therefore (all entries are to be multiplied by 

ilh) 


(■.j)l4AI 

AtUt 

Aitn 


(bj) 


0 

0 

tij 


c. Zeeman term 

The Zeeman term can be written as 


8. Evaluating tables 

a. Chemical potential 
The chemical potential can be written as 

= XI ftxJx4x(^^y (B43) 

iX ijXX' 

where 

fiXjX' = <^y^AA'Iii- (B44) 

The table is therefore (all entries are to be multiplied by 

ilh) 



4 mt 

Aim 

W 47 T 

^ 4.4 

ihi) 

hi 

0 

0 

hi 


=X ■ '^)aA' cIaC^A = X /ajA'cLffi-V, 

iX ijXX' 

(B47) 

where 

fYxjX' = <^y ^AA'I4 (n • cr)^^, . (B48) 

The table is therefore (all entries are to be multiplied by 


(•.j)l4SI 

Amt 

mil _1 

Aim 

\ ^4.4_1 

^4,7T 

(b*) 


14 (rix - iriy) 14 {nx + inx) 

-V^n^ 


d. Rashba spin-orbit interaction 
The Rashba spin-orbit interaction can be written as 


Note: Although the chemical potential gives rise to 
some non-zero terms, it does not give rise to a contribu¬ 
tion to any of the sink-source or current operators. That 
it does not contribute to the currents is clear from the 
fact that it only connects terms on the same site. To see 
that it also does not contribute to the sink-source term 
we note that the chemical potential term is rotationally 
invariant and can be considered to already be written in 
the CT-basis for any given a. That is, the A’s in the table 
can be co nsidered to be written directly in the tr-basis. 
Eq. (B20) reveals that only A A' enters 

into the sink-source term, all of which are zero in the 
table above. 


b. Kinetic term 

The spin-degenerate kinetic energy can be written as 

= X cIaCjA = X 4xoX' 4aCjA', (B45) 

ijX ijXX' 


f^SO _ Sij+x) ( 

ijXX' 

+i ~ ^i,j+y) \\i) C^Affi'A' 

= X /aiA'cLffiAS (B49) 

ijXX' 

where 

fiXjX' —^ij — Si,j+x) i~'^^y)XX' 

+* {^i,j-y ~ ^i,j-i-y) {^x)x\i) ■ (B50) 

The table is therefore (all entries are to be multiplied by 


(bj) 

aTTTT 

^4.4 

Atm 

m.it 

Aim 

^4.4 

4.4 

(i, i-\- x) 

0 


Oii,i-\-x 

0 

(i, i — x) 

0 


^i,i—x 

0 

(b* + y) 

0 



0 

(b* - y) 

0 

'^OLi^i — y 


0 
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e. Total 


that Eq. (B231 has to be modified to read 


For a total Hamiltonain of the form 

H = + H^+ (B51) 


the total table is obtanied as the sum of the individual 
tables (all entries are to be multiplied hy i/h) 


ihj) 


Atm 

m.it 

Aitn 

4 tut 

^4.7t 

(s,s) 

VzUz 

14 - iriy) 

14 + iriy) 

-VzUz 

(s — X, s) 

^s—x,s 


^s — x,s 

^s—x,s 

(s + X, s) 


C^s-\-x,s 

^S+X,S 

^S-\-X,S 

{s - y,s) 

^s—y,s 

iO-S — y^S 

y,s 

^s — y,s 

{s + y,s) 




^s-\-y,s 


We have here ignored the contributions from because 
as noted above all relevant terms are zero. We have also 
assumed that tij only is non-zero for nearest neighbour 
terms. 


^ - E [Jn + ^n) + I J . 

(B54) 


Taking the expectation value of this expression and using 
that As = -14c(cs^Cst) = (-^sc(Cs^c^;)) = (A*)* we 
find 


^ =(^*) - E ■ (B55) 


Rewriting this as 


9. Note on superconductivity 

All terms in the Hamiltonian treated so far are 
quadratic of the form cl^Cj\. However, even within a 
mean-field treatment of ordinary s-wave superconduc¬ 
tivity, we will also encounter terms in the Hamiltonian 
which are of the form 


-,SC 




,t J 
'it 4 


Cit ^it 7 


(B52) 


where A^ = —VsdciiCi^), for some 14c- This term is a 
spin-singlet and therefore appears similarly in any given 
spin-basis, such that we can assume that a ='[. If the 
term is included in the Hamiltonian, then it is clear from 


h I 


a , Cj,^Cst 


= - (^A*CsiCs^ - 




) (B53) 


d{Pst) = (S. 



(B56) 


it is possible to view any non-zero term as a possible 
generator of density fluctuations. Although the super¬ 
conducting condensate can carry current, it is clear from 
the vanishing expectation value that it does not in itself 
generate density fluctuations. In fact, it is clear from 
Eq. (|B54|) that the superconducting term is unable to 


generate any currents at all, including currents which 
preserves the density, as it gives rise to a purely on¬ 
site contribution to dps^/dt. The superconducting term 
is therefore not responsible for generating any currents; 
the condensate only carries them once they exist. The 
treatment of only ordinary quadratic terms is therefore 
sufficient for calculating the current also in the presence 
of s-wave superconductivity. 


















